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ABSTRACT 



Context. Recent observations of the X-ray emission from T Tauri stars in the Orion nebula have shown that they undergo 
frequent outbursts in their X-ray luminosity. These X-ray flares are characterised by increases in luminosity by two orders of 
magnitude, a typical duration of less than one day, and a significant hardening of the X-ray spectrum. 

Aims. It is unknown what effect these X-ray flares will have on the ionisation fraction and dead-zone structure in protoplanetary 
disks. We present the results of calculations designed to address this question. 

Methods. We have performed calculations of the ionisation fraction in a standard Q-disk model using two different chemical 
reaction networks. We include in our models ionisation due to X-rays from the central star, and calculate the time-dependent 
ionisation fraction and dead-zone structure for the inner 10 AU of a protoplanetary disk model. 

Results. We find that the disk response to X-ray flares depends on whether the plasma temperature increases during flares 
and/or whether heavy metals (such as magnesium) are present in the gas phase. Under favourable conditions the outer disk 
dead-zone can disappear altogether, and the dead-zone located between 0.5 < R < 2 AU can disappear and reappear in phase 
with the X-ray luminosity. 

Conclusions. X-ray flares can have a significant effect on the dead-zone structure in protoplanetary disks. Caution is required 
in interpreting this result as the duration of X-ray bursts is considerably shorter than the growth time of MHD turbulence due 
to the magnetorotational instability. 

Key words, accretion, accretion disks - MHD - planetary systems: protoplanetary disks - stars: pre-main sequence 



1. Introduction 

Observations of young stars in a variety of star forming 
regions, ranging from the trapezium cluster in Orion to 
the Taurus-Auriga complex, have shown that protostel- 
lar disks are ubiquitous (e.g. Beckwith & Sargeant 1996; 
O'Dell et al. 1993; Prosser et al. 1994). These disks of- 
ten show evidence for active accretion with a canonical 
mass flow rate onto the central star of ~ 10~ 8 M© yr _1 
(e.g. Sicilia-Aguilar et al. 2004), requiring a mechanism 
to transport angular momentum within the disks. So far 
only one mechanism has been shown to work: MHD tur- 
bulence generated by the magnetorotational instability 
(MRI) (Balbus & Hawley 1991; Hawley & Balbus 1991). 

Given that protostellar disks are cool and dense near 
their midplanes, there are questions about the global ap- 
plicability of the MRI to these disks, as the ionisation frac- 
tion is expected to be low (Blaes & Balbus 1994; Gammie 
1996). Indeed, nonlinear magnetohydrodynamic simula- 
tions of disks that include ohmic resistivity (Fleming, 
Stone & Hawley 2000) have indicated that for magnetic 
Reynolds numbers i?e m below a critical value -Re™*, tur- 



bulence is not sustained and the disks return to a near- 
laminar state whose internal stresses arc too small to ex- 
plain the observed mass accretion rates onto T Tauri stars. 

There have been a number of studies of the ionisation 
fraction in protostellar disks. Gammie (1996) first sug- 
gested that disks may have magnetically "active zones" 
sustained by thermal or cosmic ray ionisation, adjoining 
regions that are "dead-zones" where the ionisation frac- 
tion is too small to sustain MHD turbulence. Sano et al. 
(2000) examined this issue using a more complex chemical 
model that included dust grains. Glassgold et al. (1997) 
and Igea et al. (1999) examined the role of X-rays as a 
source of ionisation in protostellar disks, and highlighted 
doubts about whether Galactic cosmic rays could pene- 
trate into the inner regions of protostellar disks because 
of the stellar wind. Fromang, Terquem & Balbus (2002) 
examined the influence of gas phase heavy metals and 
demonstrated the potential importance of charge-transfer 
reactions, and Scmenov et al. (2004) studied disk chem- 
istry and the ionisation fraction using a complex reaction 
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network drawn from the UMIST data base. 

Recent observations of the X-ray emission from T 
Tauri stars in the Orion nebula using the Chandra obser- 
vatory (COUP - Chandra Orion Ultradeep Project) have 
shown that in addition to providing a characteristic X-ray 
luminosity at a level of Lx ~ 10 30 erg s -1 , young stars 
emit X-ray flares whose luminosity is ~ 100 times this 
value (e.g. Wolk et al. 2005; Favata et al. 2005). These 
flares typically last for less than a day, and are charac- 
terised by a sharp linear rise in luminosity, followed by an 
exponential decay. The typical recurrence time is about 
one week, and associated with the flares is a hardening of 
the X-ray spectrum indicating a rise in the plasma tem- 
perature in the stellar corona from k^T ~ 3 keV to typical 
values of ~ 7 keV. In this paper we address the question 
of what effect these X-ray flares have on the ionisation 
fraction and structure of dead-zones in protostellar disks. 

In a recent paper (Ilgner & Nelson 2006a) we com- 
pared the predictions made by a number of chemical reac- 
tion networks about the structure of dead-zones in stan- 
dard a-disk models. This study included an examina- 
tion of the reaction scheme proposed by Oppenheimer & 
Dalgarno (1974), and more complex schemes drawn from 
the UMIST data base (Le Teuff et al. 1996). In a follow- 
up paper (Ilgner & Nelson 2006b) we examined the role 
of turbulent mixing in determining the structure of dead- 
zones in a-disk models using chemical reaction networks 
drawn from (Ilgner & Nelson 2006a). In this paper we 
continue with our work on the ionisation structure within 
protoplanetary disks and examine the effect that X-ray 
flares have on dead-zones in a-disk models using reaction 
networks drawn from Ilgner & Nelson (2006a). 

In general we find that X-ray flares can have a fairly 
dramatic effect on the ionisation structure in disks, espe- 
cially if the X-ray spectrum hardens during flares and/or 
trace quantities of heavy metals (magnesium) are present 
in the gas phase. Our disk models can be divided into three 
distinct regions: an inner region where the disk is always 
active due to thermal ionisation; a central region in which 
the dead-zone formally decreases in depth substantially 
or disappears altogether during X-ray flares, but which 
returns to being a deep dead-zone in between flares; an 
outer region beyond R = 2 AU in which the dead-zone 
depth does not change in time, and which can become 
very thin or disappear altogether in the presence of heavy 
metals and an increasing plasma temperature during out- 
bursts. 

This paper is organised as follows. In section [2] we de- 
scribe our modelling procedure, including the disk model, 
the chemical models, and our method for simulating X-ray 
flares. In section [3] we present the results of our models, 
and discuss the effects of changing the model parameters. 
In section 0] we discuss our results in the context of tur- 
bulent protostellar disks, and in section [5] we summarise 
our findings. 



2. Model 

2.1. Disk model 

The underlying disk model considered is a standard a— 
disk. Details are given in Ilgner & Nelson (2006a) and 
references therein. To recap: the disk is assumed to or- 
bit a young solar mass star and undergo viscous evolu- 
tion. We use the a prescription for the viscous stress, such 
that the kinematic viscosity v — ctcj/Q, where c s is the 
sound speed and f2 is the local Keplerian angular veloc- 
ity. Heating of the disk is provided by viscous dissipation 
alone, and cooling by radiation transport in the vertical 
direction. The disk structure is obtained by solving for 
hydrostatic and thermal equilibrium. The disk model is 
completely specified by the mass accretion rate, M and 
the value of a. In this paper we consider a single disk 
model with M = 10~ 7 M yr _1 and a = 10~ 2 . The mass 
is 0.0087 M Q between 0.1 < R < 10 AU. 

2.2. Kinetic models 

We have applied two kinetic models to evolve the gas- 
phase chemistry. Both reaction networks have been de- 
scribed in Ilgner & Nelson (2006a), where the models 
were given the labels modell and model3, respectively. 
For continuity, we maintain this labelling convention in 
this paper. To recap: modell refers to the kinetic model 
of Oppenheimer & Dalgarno (1974) while model3 is linked 
to the more complex UMIST database. The underlying ki- 
netic scheme of the Oppenheimer & Dalgarno model, in- 
volves two elements, five species, and four reactions. These 
species are free electrons "e~", a representative molecule 
"m", a heavy metal atom "M", and their ionized coun- 
terparts "m + " and "M + ". We used the reference values 
of the rate coefficients given in Ilgner & Nelson (2006a). 
The UMIST kinetic model was constructed by extracting 
all species and reactions from a set of 174 species in the 
UMIST database containing the elements H, He, C, O, N, 
S, Si, Mg, Fe, resulting in 1965 reactions being included. 
A detailed description of both kinetic models is given in 
Ilgner & Nelson (2006a). 

Apart from the ionisation rate C which is discussed in a 
separate section below, all the other parameters are taken 
from Ilgner & Nelson (2006a). 

2.3. X-ray flares 

We assume that ionisation of the disk material arises be- 
cause of incident X-rays that originate in the corona of 
the central T Tauri star. We neglect contributions from 
Galactic cosmic rays as it remains uncertain whether they 
can penetrate into inner disk regions we consider due to 
the stellar wind. In our previous work (Ilgner & Nelson 
2006a, 2006b) we introduced the X-ray luminosity L x , 
which was constant. In this work the X-ray luminosity be- 
comes time dependent, L x (t) because of the X-ray flares. 



M. Ilgner, R.P. Nelson: Ionisation fraction in disks 



3 



-2 



1 2 3 4 5 6 7 

t [days] 

Fig. 1. Time sequence of the X-ray flares #C/Co- The per- 
turbation is periodic in time with a periodicity r x = 
7 days while the width of a single perturbation is t a = 
12 hr. 



We now adopt a model in which the X-ray luminosity 
maintains a base value, L x = L x , on top of which are 
superposed X-ray flares with peak luminosity L x . In or- 
der to maintain the compatibility with our previous work 
(Ilgner & Nelson 2006a, 2006b) we use L x = 10 30 erg/s. 

We approximate the flare temporal morphology by a 
sequence of outbursts that arise periodically in time. The 
shape of a single flare is modelled by using a Gaussian 
profile. We note that the symmetric profile thus obtained 
is simpler than the observed flare morphology which is 
more accurately characterised as a linear rise and expo- 
nential decay, but we believe that our model captures the 
essentials of how X-ray flares affect the ionisation struc- 
ture in the disk independently of such details. We assume 
that one flare occurs per week, and each outburst lasts 
for 24 hours. These values are similar to those observed 
during the Chandra Orion Ultradeep Project (COUP) as 
reported by Wolk et al. (2005). 

For a given total X-ray luminosity L x and plasma tem- 
perature k B T one can calculate the ionisation rate due to 
X-rays at each position in the disk model. This requires 
an integraton along each line of sight through the disk 
model to the X-ray source, and our method for this is 
described in Ilgner & Nelson (2006a). In this paper we 
consider both models in which the plasma temperature 
remains constant, and models in which the plasma tem- 
perature varies along with the X-ray luminosity. When 
the plasma temperature remains constant the attenuation 
of the X-rays is also constant, so the temporal morphol- 
ogy of the X-ray luminosity and the local ionisation rate 
are characterised by the same mathematical function. The 
perturbation <5£ in the local ionisation rate due to an X- 
ray flare is then given by the following time-dependent 
ionisation rate 



where 



sm = 



ACo exp 




if A - t a < t < A + t„ 
otherwise 



(2) 



where Co an d <5C denotes the base X-ray ionisation rate 
and its time-dependent perturbation, respectively. A is a 
factor that determines the amplitude of the perturbation, 
t a = 12 hr, A = tot x with m as an integer, and t x = 
7 days. Hence the characteristic width t a of a single flare 
is given by 



to 



(3) 



with <5C/A = 10~ 5 for \t— A| = t a . A time sequence of the 
modulation of 5C, with A = 100 is shown in figure ^ 

Observations indicate that the plasma temperature 
k^T increases as the X-ray luminosity does during a flare. 
This means that we need to take account of the fact that 
the penetration depth also varies with time, due to the 
hardening of the X-ray spectrum. As we will see in sec- 
tion 13.81 this can have quite dramatic effects on the lo- 
cal ionisation rate. To model this we consider a minimum 
plasma temperature T coo i which applies when the X-ray 
luminosity is at its base value L x and calculate the ion- 
isation rate ^(R, z) at each position in the disk. When 
the X-ray luminosity has reached its peak value L x we 
assume that the plasma temperature has reached its max- 
imum value Thot, and calculate the ionisation rate Cp(-R) z ) 
at each position in the disk. The time dependent pertur- 
bation to the local ionisation rate is then given by 



S((t) = Cp exp 




(4) 



at) = Co + 



(i) 



if \t — A| < t a , otherwise = 0. In the models with 
varying plasma temperature, we adopt kBT coo \ = 3 keV 
and fceThot = 7 keV, which are close to the observed values 
(Wolk et al. 2005; Favata et al. 2005) 

3. Results 

We assume L x = 10 30 erg/s for the base value of the 
X-ray luminosity, in basic agreement with recent obser- 
vations. During the X-ray flare the X-ray peak lumi- 
nosity i x is significantly larger then L x . We assumed 
L x = 100 x L x which is in line with the observed changes 
reported in Fatava et al. (2005). The ionisation rate C is 
then given by equations Q and (0 with A = 100 if k^.T 
is constant, or equations Q and if k^T varies. The 
time averaged value over one complete seven day cycle is 
<C> = Ccff = 5.81 x C when k^T has a constant value 
of 3 keV. 

We have evolved the disk chemistry using the kinetic 
models model 1 and model3 by taking the perturbation S( 
of the ionisation rate due to X-ray flares into account. 
As in our previous studies we wish to determine which 
parts of the disk are sufficiently ionised for the gas to 
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Fig. 2. modell - Column densities of the whole disk (solid line) and of the active zones (dashed and dotted lines) - 
corresponding to magnetic Reynolds numbers Re m greater than 100 - for two different heavy metal abundances Xm- 
The left panel is for im = and the right panels is for im = 10~ 12 . The upper and lower dashed lines refer to 
simulations with constant ionisation rate C = Co and C = 100 x Co while the long-dashed line refers to simulation with 
C = Ccff- The dotted line refers to the column density of the active zone at t — 10,000 yr obtained by assuming a 
time-dependent ionisation rate ((t) = Co + $C(t)- The time of the plot is at the mid-point between two X-ray flares. 



be well coupled to the magnetic field, and thus able to 
maintain MHD turbulence, and which regions are too neu- 
tral for such turbulence to be maintained. Again, we re- 
fer to those regions as being "active" and "dead" zones 
respectively, with the region bordering the two being the 
"transition" zone. The important discriminant that deter- 
mines whether the disk is active or dead is the magnetic 
Reynolds number, i?e m , defined by 



where H is the disk semi-thickness, c s is the sound speed, 
and \i is the magnetic diffusivity which is a function 
of the local free electron fraction. We adopt a value of 
Re" [t = 100 in this paper, following the value used in 
our previous publications and Fromang et al. (2002). We 
are able to calculate the distribution of Re m within our 
disks. Regions with Re m < 100 are deemed to be mag- 
netically dead, and those with i?e m > 100 magnetically 
active. The quantity x crit [e _ ] denotes the ionisation frac- 
tion x[e~] along the transition zone. 

Our approach to modelling the chemistry in disks with 
X-ray flares is as follows. First, we evolved models keep- 
ing the ionisation rate £ constant in time. We considered 
three cases with C = Co, C = Ccff' an< ^ C = 100 x Co 
The kinetic equations are solved for a time interval of 
At = 100, 000 yr. Hence, the ionisation fraction x[e~] is 
a function of time t, and in principle so is the location of 
the transition zone. However, for these models which as- 
sume a time-independent ionisation rate C, the change in 
the vertical location of the transition zone at all cylindri- 
cal radii in the computational domain was below the grid 
resolution for t > 10, 000 yr. 

Using the abundances obtained from the C = Co model 



at t = 100, 000 yr as initial abundances, we also cal- 
culated the ionisation fraction x[e~] by considering a 
time-dependent ionisation rate C — Co + $C The ki- 
netic equations here were solved for a time interval of 
At = 10, 000 yr. We restricted the time integration of 
the kinetic equations by assuming a finite maximum ab- 
solute step size A/i max = 1 hr in order to resolve the flares 
in time. Because of the adopted value of the perturbation 
A = 100, models with time-independent ionisation rates 
C = Co an d C = 100 x Co represent the limiting cases for 
models with C = Co + ^C- 

3.1. Snapshots at t = 10,000 yr with k^T constant 

We begin by discussing models which have a constant 
plasma temperature k^T = 3 keV but varying X-ray lu- 
minosity. The effects of allowing the plasma temperature 
to vary are discussed in section Rj. 31 The discussion below 
relates to the disk properties at a single point in time 
midway between two X-ray flares. Discussion about the 
time dependent behaviour of the ionisation fraction and 
dead-zones is presented in section 13.21 

Oppenheimer &: Dalgarno model 

The results obtained for modell are presented in fig. 
which shows the column density of the whole disk plotted 
as a function of radius using the solid line, and the 
column density of the active zone using either dashed 
lines (for which C = Co, C = C c ff an d C = 100 x Co) or 
dotted lines referring to C(*) = Co + <5C(^) - The dotted 
line corresponds to a time which is halfway between two 
X-ray flares. The left panel shows cases for which the 
heavy metal abundance xm — and the right panel shows 
cases with im = 10~ 12 . 
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Fig. 3. model3 - Column densities of the whole disk (solid line) and of the active zones (dashed and dotted lines) 
- refering to magnetic Reynolds numbers i?e m greater than 100 - for two different elemental abundances XMg- The 
left panel is for XMg = and the right panel is for xm% = 10~ 8 . The upper and the lower the dashed lines refer to 
simulations assuming a constant ionisation rate C = Co and C = 100 X Coj while the long-dashed dotted line refers to 
simulation with £ = £ e g. The dotted line shows the column density of the active zone at t = 10,000 yr assuming a 
time-dependent ionisation rate = Co + <^CW- The amplitude A of the perturbation SC/Co is A = 100 caused by 
increasing the X-ray luminosity — 100 x at the flare peak. 



As expected, figure [2] shows that a time-independent 
increase in the X-ray luminosity leads to a corresponding 
decrease in the depth of the dead-zone (note that the disk 
inner regions are fully active because of thermal ionisation 
of potassium [Ilgner & Nelson 2006a]). Interestingly, even 
an increase in X-ray luminosity by a factor of 100 is 
insufficient to fully ionise the disk, and a dead-zone 
remains beyond R > 0.5 AU. When we consider the 
X-ray flaring model, we find that the behaviour of the 
ionisation fraction and dead-zone depends on radial po- 
sition within the disk. The region between 0.5 < R < 1.2 
AU has the same dead-zone structure as the model with 
constant X-ray luminosity set at the base level L^- As 
will be discussed in more detail later on, this is because 
the recombination time in this region is shorter than the 
period of the X-ray outburst cycle, so the disk ionisation 
fraction remains in phase with the X-ray luminosity. 
Further out in the disk beyond R > 2 AU, we find that 
the dead-zone structure is very close to that obtained 
when the X-ray luminosity takes a constant value that 
is equal to the time-averged value of the X-ray flaring 
model (i.e. £ = Ccff)- This is because the recombination 
time is now longer than the period of the X-ray flaring, 
so the response of the ionisation fraction lags the instan- 
taneous ionisation rate. Over long evolution times the 
disk responds to the average ionisation rate. 

The right-hand panel in figure |3 shows the effect of 
introducing a small abundance (% = 10~ 12 ) of heavy 
metals. As described in Fromang et al. (2002) and Ilgner 
& Nelson (2006a), the introduction of heavy metals to 
the gas phase is expected to increase the free-electron 
fraction because of charge-transfer reactions with molec- 
ular ions. Figure [5] shows that an increase in constant 



X-ray luminosity by a factor of 100 causes the dead-zone 
to disappear beyond 2 AU. A constant ionisation rate 
Ceff corresponding to the time-average of the flaring rate 
leads to effective removal of the dead-zone beyond 3 AU. 
When we consider the X-ray flaring model we observe 
similar behaviour to that in the model without heavy 
metals. Interior to about 1 AU the dead-zone structure 
at the time of the snapshot is essentially the same as 
the one obtained when the ionisation rate has a constant 
value C = Co- As discussed in Ilgner & Nelson (2006a), 
the recombination of free electrons in this region remains 
dominated by molecular ions even when heavy metals 
are present, so this result is expected. Further out in the 
disk the recombination becomes dominated by the heavy 
metal ions, M + , and the recombination time is longer 
than the period of the X-ray flares. In this region the 
disk again responds to the time dependent X-ray flares 
as if the ionisation rate were equal to the time averaged 
value, resulting in the effective removal of the dead-zone 
beyond R ~ 3 AU. 

UMIST model 

The results obtained for model3 are presented in fig. |3| 
Here, the column density of the whole disk is plotted as 
a function of radius using the solid line while the dashed 
lines correspond to the column density of the active zone 
obtained for models with £ = ( , ( = Ceff, and ( = 100xC - 
In addition, the column density of the active zone for 
the model employing the time dependent ionisation rate 
£ = Co + is shown by the dotted line. The left panel 
corresponds to models in which the elemental abundance 
of magnesium XMg = 0, the right panel shows results for 
which XMg = 10~ 8 . 
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Fig. 4. modell - The change in the ionisation fraction x[e~]/x cr it[e _ ] is plotted against the relative change in the 
ionisation rate (Co + #C)/Co at different radial positions R along the transition zone assuming im = 0. ft is a limit 
cycle since the orbit is closed (the arrow points forward in time). In addition, the equlibrium values a^fe - ] are shown 
obtained for the unperturbed models with C = Co: and ( — 100 x Co- 



Our previous work that compared the results of dif- 
ferent chemical networks has already established that the 
more complex chemical models based on the UMIST data 
base predict deeper and more extensive dead-zones than 
the simpler Oppenheimer & Dalgarno models (Ilgner & 
Nelson 2006a) . The primary reason is the larger number of 
molecular ions that occur in the complex chemistry, lead- 
ing to a faster recombination of free electrons. Comparing 
figures ^ an d El confirms this. Considering the response of 
mode 13 to differing constant X-ray fluxes, we observe that 
very similar trends arise as in modell. For example, an in- 
crease in X-ray luminosity by a factor of 100 still leaves a 
substantial dead-zone beyond R > 0.5 AU. 

Considering model3 with X-ray flares, we see that the 
dead-zone between 0.5 < R < 2 AU is of similar depth 
to that generated by the base value of the X-ray lumi- 
nosity L x (with associated ionisation rate Co)- This arises 
for similar reasons to those given when describing modell 
above: the local recombination time is shorter than the 
time period between X-ray flares, ensuring that the in- 
stantaneous ionisation fraction remains more or less in 
phase with the X-ray luminosity. Conversely, the dead- 
zone depth beyond R > 2 AU approaches that predicted 
by the model whose ionisation rate Ceff is the time average 
of the X-ray flaring model. This arises because the recom- 
bination time here is longer than the time period between 
X-ray flares, ensuring that the local ionisation fraction re- 
sponds to the average X-ray luminosity rather than the 
instantaneous value. 

The right panel shows the effect of adding an elemental 
abundance of magnesium xyig = 10~ 8 . For a constant X— 
ray luminosity increased above the base value by a factor 
of 100, the dead-zone disappears beyond R > 2 AU. The 
dead-zone is predicted to be very shallow in this region 
by the X-ray flaring model. In both these cases, however, 



a significant dead-zone remains between 0.5 < R < 2 AU 
because the dominant sources of recombination here are 
molecular ions rather than magnesium ions. 

3.2. Time-dependent results beyond t = 10, 000 yr 
with k^T constant 

In the previous section 13.11 our analysis focused on the 
chemical state of the disk at a particular point in time 
midway between two X-ray flares after 10,000 yr of 
evolution. Since the ionisation rate £(t) = Co + <K(0 
changes periodically in time, so does the size of the active 
zone. We now focus on this periodic time dependence. 
We consider the time-dependent ionisation fraction at 
a few well-defined positions within the disk, along the 
transition zone defined by figures [2] and [3 We refer to the 
ionisation fraction at the transition zone (defined where 
the magnetic Reynolds number = 100) as the critical 
value, and denote it by x C rit[e - ]. Again we discuss the 
results of the Oppenheimer & Dalgarno model prior to 
the UMIST model. 

Oppenheimer &: Dalgarno model 

We found that the change in x[e~] with time t is given by 
a limit cycle. The limit cycle becomes apparent by plot- 
ting the ionisation fraction relative to the critical value 
in the transition zone x[e _ ]/x crit [e _ ] against the change 
in the ionisation rate (Co + <^C)/Co- The ionisation fraction 
repeats with a periodicity of r x resulting in a closed or- 
bit. Examples of limit cycles obtained for modell in the 
absence of heavy metals are shown in figure note that 
the arrow points forward in time such that the cycle is 
traced in an anticlockwise direction. The left hand panel 
shows the variation of x[e~] at R = 1 AU for modell. The 
limit cycles obtained at R < 1 AU with x M = 10~ 12 are 
quite similiar to those obtained with x M = since the 
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Fig. 5. model3 - The change in the ionisation fraction x[e~]/x cr it[e~] is plotted against the relative change in the 
ionisation rate (Co + ^C)/Co near the transition zone at R = 1 AU for different elemental abundances iM g . The left 
panel is for XMg = and the right panel is for iM g = 1CP 8 . It is a limit cycle since the orbit is closed (the arrow 
points forward in time). Note that cc[e - ]/x cr jt[e - ] is slightly less than unity when C — Co because the position of the 
transition zone does not coincide precisely with a grid cell in this model. 



molecular ion still dominates the recombination of free 
electrons there. At locations where metal ions dominate 
when they are included, such as at 3 AU (as shown in the 
right panel of figure^ the change in the ionisation fraction 
x[e~] across the limit cycle is tiny and would appear as 
a straight line when using the same scale of figure The 
results of model 1 at 3 AU without heavy metals shows 
a modest rise in ionisation fraction during the cycle, as 
shown by the right panel of figure 0| 

The time required for the system to establish a limit 
cycle varies depending on the local position along the tran- 
sition zone, but to within an accuracy of 0.1% limit cycles 
were achieved throughout the disk for model 1, with and 
without heavy metals. 

We differentiate between limit cycles for which the 
change in ionisation fraction x[e~] lags substantially be- 
hind the corresponding ionisation rate C(i), and those for 
which x[e~] does not. The former behaviour is illustrated 
by the right panel of figure 0] where the value of x[e~] 
continues to rise during the rise and fall of the perturbed 
ionisation rate (Co + ^C) /Co during an X-ray flare. This 
phenomenon occurs only in those regions where the recom- 
bination time is sufficiently long that the instantaneous 
recombination rate is smaller than the instantaneous per- 
turbed ionisation rate for the duration of the flare. Once 
the flare has died away then the ionisation fraction slowly 
decreases back down to its original value, just in time for 
the next flare to begin. This slow decrease in x[e~] arises 
because the recombination time is longer than the time 
period between X-ray flares. 

The left panel of figure 01 shows an example where the 
changes in ionisation fraction are more in phase with the 
perturbation to the ionisation rate. This occurs in regions 
where the recombination time becomes short compared 
to the perturbed ionisation rate during a flare. The on- 



set of an X-ray flare leads to a quite dramatic rise in 
the ionisation fraction at the transition zone. As the flare 
reaches its peak and begins to subside, the recombina- 
tion rate starts to exceed the ionisation rate and the ion- 
isation fraction decreases while the flare subsides. At the 
point where the flare ends the perturbation to x[e~] has 
dropped from a peak value of x[e~] ~ 8x cr it[e~] down to 
x[e~] a 4x cr it[e - ]. Once the flare has subsided completely 
the ionisation fraction drops to the value corresponding to 
the steady state obtained when the X-ray luminosity has 
its base value (with associated ionisation rate Co)- This 
final drop in the ionisation fraction to the critical value 
cc[e~ it ] at the transition zone occurs because the recombi- 
nation time is shorter than the time period between X-ray 
flares. 

Figure 01 also shows that the duration of the X-ray 
flares is too short to reach the ionisation fraction obtained 
for modell with a constant ionisation rate C = lOOCo- 
Indeed we found that for the unperturbed models the 
transition from state Xoo[e~} Co] to i M [e^; 100 x Co] occurs 
on time scales between 1 day (at R — 1 AU) and 90 days 
(at R = 10 AU), where Xoo denotes the steady fractional 
abundance. 

UMIST model 

We now continue by discussing the time dependent evolu- 
tion of the ionisation fraction x[e~] obtained for mode 13 at 
positions in the immediate vicinity of the transition zone 
for times beyond 10,000 yr. Once again we found limit cy- 
cle behaviour throughout the disk, but in order to obtain 
periodic orbits beyond radii R > 2 AU we had to evolve 
the disk slightly beyond 10,000 yr (since a steady-state had 
not been reached by 10,000 yr). Figure El shows examples 
of the limit cycles obtained at R = 1 AU in the absence of 
magnesium (left panel) and with an abundance of magne- 



8 



M. Ilgner, R.P. Nelson: Ionisation fraction in disks 




Fig. 6. - model 1 - Column densities of the whole disk (solid line) and of the active zones (dashed and dotted lines) 
- refering to magnetic Reynolds numbers i?e m > 100 - for im = 0. The dashed lines refer to simulations obtained by 
assuming unchanging ionisation rates C = Co f° r different plasma temperatures. 

left panel: The dotted line refers to the column density of the active zone at t = 10, 000 yr obtained by assuming a 
time-dependent ionisation rate caused by periodic variation of the plasma temperature between fcB2c 00 i = 3 keV and 
k B T hot = 7 keV. 

right panel: The dotted line refers to the column density of the active zone at t = 10, 000 yr obtained by assuming 
a time-dependent ionisation rate £(i) caused by both the X-ray luminosity and the plasma temperature increasing 
during flares. The disk becomes entirely active when the X-ray luminosity is increased by a factor of 100 and the 
plasma temperature is maintained at a constant value of k^T = 7 keV. 



sium xug = 10~ 8 (right panel). For that region, metals do 
not have a dramatic effect on the ionisation fraction x[e~] 
because molecular ions still dominate the recombination 
process of free electrons (see Ilgner & Nelson 2006a). It is 
clear that in both cases the free electron fraction remains 
close to being in phase with the changing X-ray luminos- 
ity, with x[e~] increasing as (Co + ^C)/Co increases, and de- 
creasing as (Co + ^C)/Co does. Here the recombination rate 
becomes larger than the ionisation rate immediately after 
the peak of the X-ray flare, causing the ionisation fraction 
to decrease as the flare intensity diminishes. The ionisation 
fraction returns to the value obtained for a steady X-ray 
flux with ionisation rate Co shortly after the flare finishes, 
and remains at this value until the onset of the next flare. 
As was the case with the Oppenheimer & Dalgarno model, 
this fall of the ionisation fraction to the base value occurs 
because the recombination time is shorter than the time 
period between X-ray flares. 

For regions beyond R > 2 AU we evolved the chem- 
istry for times beyond t = 10, 000 yr until limit cycle be- 
haviour was obtained. These limit cycles were very simi- 
lar to those shown in the right panel of figure 0] having 
similar amplitudes for the case without magnesium, and 
very small amplitudes when magnesium was included with 

JMg = 10~ 8 . 

3.3. Results with varying plasma temperature ksT 

In the preceding sections we discussed the effects of chang- 
ing the X-ray luminosity during a flare on the ionisation 



structure of a disk model, while keeping the plasma tem- 
perature constant. Observations, however, indicate that 
the plasma temperature ksT increases during flaring ac- 
tivity, resulting in the hardening of the X-ray spectrum 
and an increased ability of X-ray photons to penetrate 
into the disk (Wolk et al. 2005). Modelling of X-ray 
sources often requires a two-phase emission model (Mewe, 
Kaastra & Liedahl 1995), corresponding to two different 
temperatures for the cooler and hotter plasma components 
of the corona. Observations indicate that the cooler com- 
ponent remains essentially unaffected during a flare, while 
the hotter component increases significantly. In our model 
we assume that fcelcooi = 3 keV and fceThot = 7 keV, in 
close agreement with observations. 

The requirement that we temporally resolve the X-ray 
flares while solving the kinetic equations for the chem- 
istry means that a maximum time-step of 1 hour was 
adopted. This causes the calculations to become very ex- 
pensive computationally, such that for the UMIST model 
model3 a calculation lasting for 10,000 yr requires a run 
time of approximately three weeks for each grid-point of 
the disk model. As a consequence we have been forced to 
only consider the Oppenheimer & Dalgarno model mode 11 
in the following sections. We simply note at this point 
that model 1 tends to generate thinner dead-zones than 
model3, such that the results described below represent 
an optimistic view of which fraction of the disk can sus- 
tain MHD turbulence. We note further, however, that 
similar trends in the results are obtained with model 1 
and model3 when varying physical parameters, such that 
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model 1 gives a reasonable picture of how the ionisation 
fraction obtained using the more complex model3 re- 
sponds to changes in ionisation rate etc. 

We begin by presenting the results of a model in 
which we consider the effects of increasing the plasma 
temperature during an X-ray flare, keeping the X-ray 
luminosity Lx constant. This allows us to isolate the 
effects of the hardening of the X-ray spectrum due to 
the rising plasma temperature. We then consider a model 
for which both the plasma temperature and the X-ray 
luminosity increase during flaring activity. 

Snapshot at 10,000 yr with Lx constant and fc^T 
varying 

The results for this calculation are shown in the left panel 
of figure El which shows the column density of the whole 
disk plotted as a function of radius using the solid line. 
The dashed lines show the column density of the active 
zone for two calculations that assumed constant plasma 
temperatures of UbT = 3 keV and /cgT = 7 keV, re- 
spectively, and the dotted line corresponds to the model 
whose plasma temperature rose to /c^T = 7 keV during 
flares only without an accompanying increase in Lx above 
10 30 erg s _1 . Note that these figures refer to modell with 
no heavy metal included. 

Figure [7| shows contours of the relative ionisation rate 
when the plasma temperature is raised from ksT = 3 to 
7 keV. The change in £ is dramatic in the deeper, more 
shielded regions simply because of increased penetration 
induced by the hardening of the X-ray spectrum. Near the 
disk surface the change is only slight, but in the interior 
the local ionisation rate can increase by more than a factor 
of 1000. 

The upper dashed line in the left panel of figure El 
shows that having a constantly increased value of 
ksT = 7 keV significantly reduces the depth of the dead- 
zone throughout the disk. The dotted line, however, shows 
that at a point in time midway between two peaks in the 
plasma temperature, the dead-zone is unaffected interior 
to R < 1.2 AU because of the fast recombination time 
there. Beyond R > 2 AU the dead-zone is a little thinner 
where the recombination time is longer than the period 
between increases in plasma temperature, corresponding 
closely to that which would be obtained by exposure to the 
time averaged value of the time dependent ionisation rate. 

Results with Lx and fc^T varying 
A snapshot of the results obtained after 10,000 yr when 
the X-ray luminosity and the plasma temperature increase 
during a flare is shown in the right panel of figure El The 
solid line gives the column density of the whole disk, and 
the dashed line gives the column density of the active zone 
obtained when the X-ray luminosity is constant and takes 
its base value L x . The disk becomes active everywhere 
if we take a constant ionisation rate with Lx = 100L X 
and k-gT = 7 keV. The dotted line corresponds to the 
model for which the plasma temperature and X-ray lumi- 
nosity increase during a flare. Note that these models are 
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Fig. 7. The ratio of the effective X-ray ionisation rates 
Co, 7 kcv/Co,3 kcv- Co, 3 kcv and Co, 7 keV refer to local ion- 
isation rates for plasma temperatures k B T = 3 keV 
and k B T — 7 keV, respectively. The disk parameter are 
a = 10~ 2 and M = 10~ 7 M Q yr _1 . The contour lines refer 
to values of Co,7 kcv/Co,3 kcv: 10, 10 2 , and 10 3 . 



modell in the absence of heavy metals, and the snapshot is 
taken at a time that is midway between two X-ray flares. 
In the inner regions R < 1.2 AU the dead-zone size is 
again unaffected by the flares between outbursts because 
of the rapid recombination time there. Beyond R > 2 AU, 
however, the dead-zone depth has been decreased sub- 
stantially such that the column density of the active zone 
includes approximately 80 % of the matter in the region 
between R = 3 - 7 AU. Although not plotted in figure El 
the addition of heavy metals with im = 10 -12 causes the 
dead zone to disappear completely beyond R > 2 AU, al- 
though a dead-zone ranging between 0.4 < R < 2 AU 
remains. 

We now consider the time dependent evolution after 
10,000 yr has elapsed. The evolution of the ionisation frac- 
tion at the disk midplane (and not at the transition zone 
as in previous similar figures) as a function of the per- 
turbed ionisation rate is shown in figure |S1 at radii R = 1 
AU and R = 5 AU for modell in the absence of heavy met- 
als. The first thing to note is the large perturbation to the 
ionisation rate experienced at the disk midplane, as illus- 
trated by the range of values for (Co + <$C)/Co- Second, the 
ionisation fraction, relative to the critical value x C rit[e _ ], 
at the midplane located at R = 1 rises from 0.1 (imply- 
ing the region is dead) to a value of sa 2 (implying that 
the region is active). Thus, the inner regions of the disk 
show the limit cycle behaviour, and oscillate from having 
a dead-zone near the midplane to having an active zone 
there which lasts for a duration slightly longer than half 
the flare (i.e. about 12 hours). This region remains for- 
mally "dead" for the other six and a half days during the 
X-ray flaring cycle. 

Although at first sight this would appear to indicate 
that the stress associated with MHD turbulence in a pro- 



10 



M. Ilgner, R.P. Nelson: Ionisation fraction in disks 



R = 5 AU 
c» = 



100 

1 + <5<r/<r= 



10000 



100 10000 

1 + &</<<. 



Fig. 8. modell - The change in the ionisation fraction x[e _ ]/x cr it[e _ ] is plotted against the relative change in the 
ionisation rate (Co + <5C) / Co at different radial positions R with z/R = (midplane). x cr it[e _ ] denotes the ionisation 
fraction at the transition layer. Each graph describes a limit cycle since the orbit is closed (the arrow points forward in 
time). The perturbation <5C/Co m the ionisation rate is based on the assumption that both the X-ray peak luminosity 
Lx = 100 x -^x an d the temperature of the hot plasma component k B T hot — 7 keV at the flare peaks. 



toplanetary disk would be episodic, leading to episodic 
mass accretion, it is important to bear in mind that the 
local growth rate for the MRI is approximately the local 
dynamical time (Balbus & Hawley 1991), which is longer 
than 12 hours in the regions beyond 0.5 AU. It is there- 
fore unclear whether such a region of the disk could main- 
tain fully developed turbulence and an associated dynamo 
when perturbed by X-ray flares of short duration. This is- 
sue is discussed below. 

The right panel of figure |S1 shows the evolution of 
a;[e _ ]/a; cr it[e _ ] versus the perturbation to the local ioni- 
sation rate (Co + <5C) / Co at the disk midplane at R = 5 
AU. Here it is clear that limit cycle behaviour is again 
obtained, with the ionisation fraction showing a closed 
orbit. In this part of the disk, there remains a thin dead- 
zone throughout the cycle when heavy metals are not 
present. The addition of heavy metals with abundance 
%Mg = 10 -12 renders the disk fully active in this region 
for the whole duration of the X-ray cycle when both the 
luminosity and plasma temperature vary. 



4. Discussion 

In this paper we have calculated the response of the ioni- 
sation fraction in a standard a-disk model to the presence 
of X-ray flares originating in the corona of the central T 
Tauri star. For comparison purposes we have examined 
how the disk responds to constantly acting X-ray induced 
ionisation rates Co (the base ionisation rate), Ceff (the time 
averaged ionisation rate), and 100Co (which is the peak 
ionisation rate during a flare). We have considered cases 
where the plasma temperature retains a constant value of 
k-gT — 3 keV during the flare cycle, and cases where the 
plasma temperature increases from 3 keV to 7 keV at the 
peak of a flare. Broadly speaking, the behaviour of the 



dead-zones in our disk model can be divided into three 
basic regions when subject to X-ray flares: 
(«). An inner region R < 0.5 AU where the disk is entirely 
active due to thermal ionisation of potassium throughout 
the calculations. 

(ii). A central region 0.5 < R < 2 AU where the disk 
exhibits a significant dead-zone between outbursts, but 
whose dead-zone can be significantly diminished or be 
removed altogether near the peak of an outburst. This 
behaviour occurs because of the short recombination time 
here, even when heavy metals are included in the gas phase 
due to the recombination being dominated by molecular 
ions. This causes the ionisation fraction to remain close 
to being in-phase with the ionisation rate during and be- 
tween flares. Our calculations indicate that allowing the 
plasma temperature to increase to k^T — 7 keV during 
X-ray flares can cause the magnetic Reynolds number to 
exceed 100 for about 12 hours after the peak of the flare, 
such that the dead-zone here is formally removed for that 
time. 

(Hi). An outer region R > 2 AU where the dead-zone size 
and depth does not change significantly during the flar- 
ing cycle. Here the dead-zone size corresponds closely to 
that obtained by the time-averaged X-ray flux, because 
the recombination time exceeds the time period of the 
X-ray outburst cycle. The depth of the dead-zone size 
in this region is significantly diminshed if heavy metals 
are present in the gas phase or the plasma temperature 
increases during flares. The dead-zone may disappear al- 
together if metals are present and the plasma temperature 
increases during flares. In the absence of heavy metals and 
an increasing plasma temperature, the time averaged ion- 
isation rate during a flare cycle rj 5.8Co, and a significant 
dead-zone remains. 



M. Ilgner, R.P. Nelson: Ionisation fraction in disks 



11 



Given this disk structure consisting of an inner region 
without a dead-zone, a central region whose dead-zone 
can change significantly or be removed during parts of 
the X-ray flaring cycle, and an outer region which may 
maintain a thin dead-zone of almost constant depth, it is 
interesting to speculate how MHD turbulence in the sys- 
tem will evolve in time. Clearly the innermost regions can 
sustain MHD turbulence continuously, as can the magnet- 
ically active component of the outer disk. The behaviour 
of the central region whose magnetic Reynolds number 
varies significantly throughout the X-ray flaring cycle is 
less clear. In this region the dead-zone remains deep for 
six and a half days and is formally reduced in size or re- 
moved for only about 12 hours. 

In the ideal MHD limit the fastest growing mode of the 
MRI has wave number /c max defined by k m3 ,-xVA ~ where 
va is the Alfven speed and f2 is the local Keplerian angular 
velocity. The growth rate associated with this mode is ~ SI 
(actually closer to 0.75^ [Balbus & Hawley 1991]). Ohmic 
resistivity strongly affects the growth of linear modes with 
wave number k when the associated diffusion rate k 2 rj ex- 
ceeds the linear growth rate. Comparing these rates for the 
fastest growing mode provides a condition for the effects 
of resistivity to dominate over the growth of the mode due 
to the MRI: 

^max 7 ? k ma:x VA (6) 

which yields the condition 

I < < 7 > 

Equation Q is often used to provide an alternative defi- 
nition of the magnetic Reynolds number Re' m — v 2 A /(riVl) 
(e.g. Sano, Inutsuka & Miyama 1998) which is related to 
the definition used in this paper by 

Rem = (l Re > m , ( 8 ) 

where (5 = 2c 2 Jv 2 a = P eas /P mag . P gas and P mag are the 
gas and magnetic pressure. Thus a value of Re m — 100 
corresponds to a value of Re' m = 1 when (3 — 200, such 
that resistivity dominates over the growth of the fastest 
growing mode for values of i?e m and Re' m smaller than 
these. 

For the purpose of illustration, let us consider a sit- 
uation where f3 = 200 in the central regions of our disk 
near the location R = 0.5 AU, where the inverse of the 
Keplerian angular velocity SI -1 = 20.5 days. An X-ray 
flare arises, raising the value of Re m > 100 for 12 hours, 
allowing exponential growth of the fastest growing mode 
for this time such that an amplification of w 2.5% oc- 
curs. Given our definition of Re m it is clear that any field 
amplification that occurred while Re m > 100 is diffused 
away within the next 12 hours after Re m drops below 100, 
and the growth of this mode is unable to amplify the field 
or drive the disk toward a turbulent state. The situation 
can in principle be different for longer wavelength modes. 



Consider the longest wavelength mode that can fit within 
the disk vertical extent with wavelength A = 2H, and as- 
sociated wavenumber fc m i n = tt/H. The wave number as- 
sociated with the fastest growing mode has fc max = 10/ H 
if {3 = 200, such that A: max /A: m i n « 3. The diffusion time 
associated with the longest wavelength mode is an order 
of magnitude longer than for the fastest growing mode 
(whereas the growth rate is approximately half of the max- 
imum value). Thus any field amplification that occurs dur- 
ing the 12 hours when Re m > 100 is diffused away over the 
next w 5 days. This suggests that an X-ray flaring cycle 
with a periodicity of less than five days could in principle 
lead to gradual field amplification over successive cycles. 
The question of whether fully developed turbulence can 
arise in such a scenario is unclear, and can only be ad- 
dressed by means of non linear simulations that explicitly 
account for periodic rises in the ionisation rates due to 
X-ray flaring and electron recombination. 

The calculations presented in this paper consider only 
gas-phase chemistry, and ignore the effects of dust. As 
such they are relevant to a stage in protoplanetary disk 
evolution when substantial grain growth has occurred and 
a dense dust layer has settled near the midplane. It is well 
known, however, that small dust grains are able to sweep 
up free electrons and substantially reduce the ionisation 
fraction (e.g Sano et al. 2000; Ilgner & Nelson 2006a), 
and the assumption of gas-phase chemistry is only really 
valid when the abundance of dust grains has been depleted 
by a factor of between 10 -4 - 10 -8 below the canonical 
concentration of 10~ 12 (Ilgner & Nelson 2006a). An open 
question that we have not addressed is what happens to 
species that are adsorbed onto the surfaces of small grains 
during X-ray flares, and in particular what happens to 
the adsorbed electrons. Najita et al. (2001) considered the 
effect of nonthermal desorption of grain mantles due to 
X-rays, but did not include the desorption of electrons 
into the gas phase. We speculate that if grain mantles 
and electrons are desorbed into the gas phase during X- 
ray flares, and are adsorbed back onto the grains during 
the time between flares, then the behaviour of the dead- 
zones throughout the disk will be similar to that already 
observed in the central regions between 0.5 < R < 2 AU. 
This is because electrons adsorb onto grains very rapidly 
throughout the disk, such that the time scale of which 
electrons will be removed from the gas phase is shorter 
than the time between flares. An analysis of this issue will 
be presented in a future publication. 

5. Summary 

We have considered the effect of X-ray flares on the ion- 
isation fraction and dead-zone structure in protoplane- 
tary disks. These flares can have a significant effect on 
dead-zones if much of the submicron sized dust in the disk 
has undergone grain growth and settled toward the mid- 
plane, if the plasma temperature increases significantly 
during X-ray flares, and/or trace quantities of heavy met- 
als (magnesium) are present in the gas phase. 



12 



M. Ilgner, R.P. Nelson: Ionisation fraction in disks 



Questions remain, however, about the disk response in 
regions where the dead-zone is removed during an X-ray 
flare, but reappears during the low state. In particular it 
is not clear that MHD turbulence can be generated there 
as the growth time of the MRI is longer than the duration 
of the flares. This issue needs to be addressed using non 
linear simulations that take account of the time dependent 
ionisation rates induced by X-ray flares. 
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